


Institutional Archive of the Naval Postgraduate School 





Calhoun: The NPS Institutional Archive 


DSpace Repository 


Theses and Dissertations 


1975 


l. Thesis and Dissertation Collection, all items 


Analysis of a high resolution deep ocean 
acoustic navigation system. 


Durham, James Leighton 


Massachusetts Institute of Technology 


http://hdl.handle.net/10945/20728 


Downloaded from NPS Archive: Calhoun 


DUDLEY 
KNOX 
LIBRARY 





http://www.nps.edu/library 


Calhoun is the Naval Postgraduate School's public access digital repository for 

research materials and institutional publications created by the NPS community. 

Calhoun is named for Professor of Mathematics Guy K. Calhoun, NPS's first 
appointed — and published — scholarty author. 


Dudley Knox Library / Naval Postgraduate School 
411 Dyer Road / 1 University Circle 
Monterey, California USA 93943 


ANALYSIS OF A HIGH RESOLUTION DEEP OCEAN 
ACOUSTIC NAVIGATION SYSTEM 


James Leighton Durham 











DUDLEY KNOX LIBRARY 
NA TGRADUATE SCHOOL 
MONTEREY, CALIFORNIA 93940 


p 


è 






- 
^ 





» 


ANALYSIS OF A HIGH RESOLUTION DEEP OCEAN 
ACOUSTIC NAVIGATION SYSTEM 


by 
James Leighton Durham 


B.S. U.S. Naval Academy 
(1970) 


Submitted In Partial Fulfillment of the 
Requirements for the Degree of 


Ocean Engineer 

at the 

Woods Hole Oceanographic Institution 
and the 

Massachusetts Institute of Technology 

and for the degree of 

Master of Science in Ocean Engineering 
at the 

Massachusetts Institute of Technology 


January, 1975 





ANALYSIS OF A HIGH RESOLUTION DEEP OCEAN 
ACOUSTIC NAVIGATION SYSTEM 


by 


James Leighton Durham 
Lieutenant, United States Navy 


Submitted to the Department of Ocean Engineering in 
partial fulfillment of the requirements for the degrees of 
Ocean Engineer and Master of Science in Ocean Engineering 


АВ5ТКАСТ 


A high resolution acoustic navigation system for ocean use 
is being developed at the Woods Hole Oceanographic Institution. 
The system can yield navigation fixes with respect to a bottom 
moored reference net with accuracies (on a fix to fix basis) of 
a few centimeters. In order to use the system to best advan- 
tage a survey 1S required to determine precisely the relative 
positions of the net elements. Each element combines a pulse 
transponder with a continuous wave (CW) beacon. Accumulated 
phase (Doppler shift of the CW beacon) between survey points is 
measured as well as acoustic travel times between survey points 
and transponders. Non-linear regression techniques are employed 
to develop a maximum likelihood estimator for net element 
positions based on these phase and travel time measurements. 

An approximate error covariance matrix is generated and an 
optimum choice of survey points is indicated. The combined 
System, using these selected locations for performing the survey, 
can yield reference mooring coordinates with error of +l meter. 
Improved precision appears to be limited by inaccuracies in the 
pulse and Doppler measuring system. | 
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I INTRODUCTION 

There are many methods currently used to position ships, 
submersibles, buoys or submerged instruments in the deep ocean. 
Broadly speaking they can be divided into two types. The first 
method positions via electromagnetic transmission from a shore 
station or satellite whose position or orbit is precisely known. 
Such systems operate at long ranges but are limited in accuracy 
to about 100-200 meters at best. Furthermore, these systems 
require above surface antennas and are therefore of little use 
in positioning entirely submerged vehicles. The second method 
positions via acoustic transmission from a set of transponders 
or beacons, usually moored to the ocean bottom, whose relative 
positions are precisely known. When operated in the pulsed or 
E ponder mode, position inaccuracies of an acoustic system 
can be as small as 10-20 meters in 5 km water depths. When 
operated in a continuous, or Doppler mode, the accuracy can 
ось 3-4 centimeters. The range of acoustic systems is 
usually severely limited by a number of factors, the most im- 
portant of which are acoustic refraction and attenuation of 
Sound in sea water due to spreading and absorption in the ocean. 
Current acoustic systems have ranges of about 10 km. 

A navigation system utilizing acoustic beacons (Fig. 1) 
which operate in both a pulse and continuous SUA mode is being 


developed. This system capitalizes on the attributes of both 
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the pulse and Doppler modes and will be capable of positioning 
a platform with respect to an array of reference beacons to an 
accuracy of 1-2 meters and of repositioning a platform within 
10 centimeters of a previous position. To make full use of the 
system capabilities a survey is required to accurately deter- 
mine the relative positions of the reference beacons. The pur- 
pose of this investigation is to develop an accurate survey 
technique. Algorithms for estimating the beacon positions and 
the errors in those poSitions are derived. Performance is 
predicted using computer models. The estimation method used is 
a straightforward extension of the pulse system technique des- 


E ed By Hunt et al. (1974). 


IJ Description of System 


Pulse (transponder) navigation systems measure the travel 
time of an acoustic pulse to estimate the slant range to a 
transponder. Shipboard processing equipment interrogates the 
bottom moored transponders and measures the time of arrival of 
SeemceP!y pulse from each transponder. Sound velocity correc- 
tions and corrections for fixed system delay times are applied, 
and travel times are converted into estimates of slant range to 
Seemetansponders. The accuracy of the measurement is depen- 


dent upon the pulse bandwidth and the signal-to-noise ratio in 


the received pulse. Systems of this type are widely in use in 
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oceanographic research, submerged rescue operations and commer- 
cial enterprises in the ocean (Boegeman et al., 1972; Cestone 
Est. George, 1972; Van Calcar, 1969). 

A continuous wave (CW) beacon system has been developed 
at the Woods Hole Oceanographic Institution to accurately 
track the position of ship-mounted and ‘submerged hydrophones 
(Porter et al., 1973). The system measures the frequency 
change of a continuous tone due to hydrophone motion.  Hydro- 
phone velocity and displacement are derived from the beacon 


Signal Doppler shift which is given by 


Af =f, - fp = с, 


where Ep 1S the Doppler shifted beacon frequency due to a 
hydrophone velocity v, Еһ is the beacon frequency and c is the 
ШОШО sound in water. Doppler shift is proportional to 
beacon frequency. Hence the accuracy of the measurement is 
dependent on the beacon frequency and the signal-to-noise 
ratio in the narrow band that encompasses the maximum Doppler 
Shift. Because the beacon frequency is much greater than the 
pulse bandwidth of transponder systems the Doppler system can 
Eure range changes much more accurately than the pulse 


system assuming that both systems have the same signal-to- 


ША Таро: In the present implementation of the system, 
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accumulated phase change ff(t) = 2 т(АЕ)Е 15 measured by digital 
techniques and hydrophone displacements of a quarter wavelength 
can be resolved. A wavelength for 13 kHz (a typical beacon 
frequency) is А = 11.5 centimeters, so that resolution is 
about А /4 23 cm. Phase changes due to inhomogeneities and 
fluctuations in the water column are negligible over the range 


of the system. We discuss this in Section III. 2. 


2) Previous Survey Techniques 


Several investigations have been conducted to develop an 
optimum technique for determining the coordinates of ocean 
bottom acoustic transponders. No satisfactory survey has been 
developed previously for an array of bottom-moored CW beacons. 

ENEUCvS of acoustic transponders generally fall into two. 
categories: baseline crossing and iterative techniques. The 
baseline crossing technique, described by Haehnle (1967) and 
mt (1967), is often referred to as the conventional or class- 
ical transponder survey method. The method requires that the 
depth of the transponders be previously determined from inde- 
pendent measurements. The baseline length, i.e., the horizon- 
tal distance between two transponders, is determined by steaming 
Peross a baseline and measuring the slant range to the two 


transponders during the traverse. This requires accurate ship 


Se 





positioning and large amounts of ship time to yield accurate 
results. More recent iterative techniques are based on or are 
КИШ ат to that described by Vanderkulk (1961). Vanderkulk 
showed that for three bottom transponders and six co-planar 
survey positions whose depths are known, the positions and 
depths of the transponders can be found by solving six linear 
equations. A unique solution exists when the survey positions 
do not lie on a conic section. When more than six survey posi- 
tions are used the additional equations result in an over- 
determined situation for which minimum least squared error 
techniques have been successfully applied, (Lowenstein, 1965; 
ШЕСІ, 1974; Mourad et al., 1972; Heckman and Abbott, 1973; 


fet al., 1974). 


3) Pulse-Doppler Survey Technique 

The pulse-Doppler survey consists of making measurements 
of travel times of acoustic pulses and accumulated phase 
changes of CW signals at various survey points. These are con- 
verted to slant ranges and slant range differences, respect- 
ively. The pulse-Doppler survey can yield more accurate results 
than the pulse survey alone due to the higher resolution of the 
Doppler system measurement. 

The slant range (SR) between a transponder and a survey 


point is equal to the travel time (TT) multiplied by an 
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effective sound speed, Ve: 

sh = Vers TT 
The effective sound speed is a function of the sound velocity 
Ese, the depths of the transponder, and the survey point 
and varies with the amount of acoustic refraction (Vass, 1966). 
ии Survey points there are m measurements of slant range 
for each transponder. 

The change in slant range, DSR, between a beacon and two 
survey points is proportional to the change in accumulated 
phase, (t5) - Ø (t) and the wavelength of the signal, 1 = c/fp: 

DE E 721 
The received signal, 

A E) expr. = 72) 1 Е-Е } 
TE f = fp +fpv/c ЕВЕ ОВР Feral t it ccm imequcncy can 
be written as 

E EZ m ft ӘЙЕЛІ! 


where the phase variation due to a finite receiver velocity is 


ШІ! - 2т Ent v/c (Porter et al., 1973). In general v isa 
function of time, v = v(t). However, we assume that the 
sampling interval is short enough to consider v constant. In 


practice, an interval of .1 sec is used so the assumption is 


valid. Measurement of accumulated phase change provides 
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velocity and displacement information along the beacon-survey 
Eu vector, 1.e., the displacement is the change in slant 
range. For m survey points there are м ~ | measurements of 
slant range differences for each beacon. 

The measurements of travel times and phase changes may 
be made at the same or different Survey points. Only the case 
where the meaSurements are made at the same points will be 
discussed in detail since this technique is more efficient in 
practice. The method of solution for the transponder coordin- 
ates is the same in either case but the associated errors will 
be different. 

Non-linear regression methods as outlined by Draper and 
Smith (1966) are used to analyze this survey problem. Under 
the assumptions given in the following section the maximum 
likelihood estimator for the beacon positions is found. An 
approximate covariance matrix for the estimates and a pro- 
cedure for optimizing the choice of survey points is developed. 
The sources of error and their effect on the estimates of the 


lmeacon coordinates are discussed. 
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II PARAMETER ESTIMATION 

tine problem of parameter estimation can be briefly des- 
cribed as follows: given a system from which information can 
be received and whose state is characterized by a set of p 
parameters, make some estimate of the state of the system from 
the information received. The information, which will gener- 
ally be perturbed by noise, may be measurements, a message, or 
some other observations on the system. 


The following definitions are used: 


Ө' = (91, 92,...,9,1 = true parameter values 
ê = OP O.» Я „ô ] = estimated parameter values 
Fi = F'(8) = [f,,f,.-..f, ] = noiseless message 
Y' = ІУІ,У2,-- Уы! = message (observations perturbed 
by noise) 
Mere o, e ШО Y denote column vectors and (') denotes 


transpose (Manasse, 1960). 

When a noiseless message F(g)is received, the information 
is sufficient to determine 0 exactly if the number of measure- 
ments n is not less than the number of parameters p. When 
n = p the problem is called the MINIMUM DATA case. When n >p 


mers Called the REDUNDANT DATA case. 
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ID The Maximum Likelihood Estimate 


When the observations Y are received, some criteria must 
be used to determine an optimum estimate 6 since š5 dQiffEeES 
E Mie to the presence of noise. One often used decision 
rule is based on the maximum likelihood principle. This 
principle prescribes that the observer choose the Ə which 
renders the observations Y most likely; that is, it maximizes 
the conditional probability density P(Y | 0). When expressed 
ШИЕ nction of Y and 8 this is called the likelihood function 
EN Sometimes written L(Y;9). When the a priori probability 
EN P(8) can be considered to be a constant over the 
region of interest then the maximum likelihood estimator also 
maximizes P(Olv). This can be shown by the use of Bayes 
т со express the a posteriori probability density 


DAS) in terms of P(Y|9) and P(P): 


P(y|@) P(9) 
P(8|Y) = - C,P(Y| 9) (1) 
P) 
where P(y) = fP(v|g) P(9) ae 


EM that P(Y) is not a function of © and can be treated as a 
Benstant of proportionality). 
tie estimate 8 which meamea he likelihood function 


І(Ү;Ө) must satisfy the equations 
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Small variance is a desired property. An asymptotic measure 
of this property is efficiency. An estimator is said to be 
efficient if the asymptotic variance O estimate is no 
larger than the asymptotic variance obtained using any other 
estimator. The maximum likelihood estimator, when it is con- 
tinuous and the first and second derivatives exist and are 
absolutely integrable, is efficient. It can be shown that 
maximum likelihood estimates are asymptotically normally dis- 
tributed with mean 9 and variance 1/R* (0) where 

9109 L : 


&^ (8)- E — (kendall and Stuart, 1967), 
98 


In the case of additive Gaussian noise the form of the 


likelihood function can be determined and expressions for the 


variance of the estimates can be developed. 


Sy the Maximum Likelihood Estimate in the Case of Additive 


Gaussian Noise 
The message is assumed to be perturbed by addıtive 
Gaussian noise, e, so that 
Y = F(8) +e (2) 
Where e is characterized by an n x n moment matrix M; 


A (ye pu) = E(ee') 


djo- 





un nHltidoxmensional Gaussian distribution of the error сап 


be expressed as 


or 

P(Y - E) 7 C2 exp{-1/2 Ix - z | y - Е | 
where Co 1s a eic constant (see, for example, p. 151 
Davenport and Root, 1958). 

From equation (2) it can be seen that the probability 
density of Y given F is equal to the probability density of 
the error; that is 

ШИН E ETE) 
Since F is a function of © this density is also the proba- 
bility density of Y given 8. The likelihood function, there- 
fore, is equal to the EIS density function of the 


error, namely a multidimensional normal distribution, that 


15, 


Е 
|< 
D 
И 


Р(У[Ө) = Р(Ү - Е(Ө)), 


pnus 


|| 


c2 exp{ 1/2[% - Е] Ми -Е]]. o 


The likelihood function is a maximum when the magnitude of 


NO) 


the exponent, ВА EJ ml [y = EM , is a minimum. 
rhe noise e has zero mean and the error in each 


measurement is statistically independent of the errors in the 


Еа 





other measurements the moment matrix M is diagonal 


J 
2 
М = (ее) = | 1 05° О 
О C2 
2 
Dr. 2 O0 
Voz 
т > 
2]: : 
Men =| o 1/0; 
L 
ore | | ch 
where denotes the variance of the error in the k measure- 


k 


ment. The moment matrix, M, can be written in terms of weight- 


ing factors w; and the standard error O . The matrix of 


-1 
weighting factors is denoted by V . 


7 
2 ^ О 2 
Maes = a I 
m | 
O Wn 
L Ж 
апа W 
1 
EC wz О 
V = O E (4) 


where м. = бб ү: 
г. 


The maximum likelihood estimate in this case is the weighted 
least squares estimate since the magnitude of the exponent 
in equation (3) is proportional to the weighted sum of squared 


Errors: 


E 





n 


Y 


= 208 i=1 (Yi 7 #1 (Ө))^ м. 


Ae 


ei rO uly - Elo) 
When the variances of the errors are equal then the weighting 
БШСЕОК5, Wj, are equal to one, and the maximum likelihood 
estimate is the familiar least squares estimate which minimizes 


E Sum of squared errors SS(B): 


2 


ss(9) = (У; - £,(9)) 


Ms 


Assume that a rather accurate estimate of 9 has been 
obtained and that F(8) can be approximated by a first order 
Taylor series 


F(9) = F(9) + х ЛӨ, 





QE | 

where X = Ен о Gorivatives with ijth 
Je 19 

Sent df; and the errors /\6 = ке б аге assumed small. 
дө; 


When this expression for F(O) is substituted in equation (3), 
amd when the a priori probability density P(@) can be consider- 
ed to be a constant over the region of interest, the probability 
density of O given Y can be shown to be a multidimensional 


Hermal distribution with moment matrix Es 


2 Ae: [7A 


р(ө|у) = Су ехр([- m. o (5) 


-21- 





where Lr -хміх- 22 осу 3), | (6) 
The moment matrix Û is known as the variance-covariance 
Matrix Since it consists of the variances (diagonal terms) and 
covariances (off-diagonal terms) of the estimates. It is also 
known aS Fisher's information matrix (Van Trees, 1968). Specif- 


ically, under the conditions of uncorrelated, zero mean Gaussian 


noise, the variance of the estimates can be expressed as 


s 2 --1 
2 2 E à f. 
Ex Ef. eo XL BI. 


: ms iE 
This estimate is approximately equal to the asymptotic 


variance 1/R* (8) о А. 


EN »Pcation of Parameter Estimation Techniques to the Survey 


Problem 

The pulse-Doppler navigation system makes measurements 
which are characterized by a set of unknown parameters. The 
unknown parameters are the geometrical coordinates of the 
beacons and the ship positions at which the measurements are 
taken. The measurements consist of slant ranges and slant 
range differences. 

An orthogonal XYZ coordinate system will be used in which 
the beacon coordinates are (0,0,2,),. (X5,0,49), and (X5, Y4, 24) 
as shown in Figure 2. (Only the 3 beacon survey problem will 


be considered here; the method can easily be extended to more 
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than 3 beacons). The depth of the Ship hydrophone is assumed 
to be a Known constant, ZS. Simultaneous measurements of slant 
ranges and slant range differences are taken at m survey points 
(Х51,Ү5., ZS) е 

The unknown parameters, that is, the beacon and ship posi- 
coo rdinates, can be expressed as a column vector 9, with 
transpose 
2t e e s 
2ر‎ 


Ө' - [0,,0 


xb X.Y XSy » YS --XS_,YS_ | 


3 Sy 1% 


where p = 2m + 6, 
The "noiseless measurements" are the exact geometrical slant 


ranges and slant range differences. These can be expressed as 


E _ b m m» ОТТО 
SR; - = SR; j (8) = eS Kad + Я Үз) + (75 Z4) ] 


where X4 - Үү = Y, = 0 i = 1,2,...,п 


5 = 1,2,3 


Similarly, the slant range differences between two survey points 


meena 11 and the aa beacon can be expressed as 


2 2 2 1072 
) + (YS; 4-¥.) + (28-24) ] 
2 
Е + (28-2. ) | 272 


DSRij (9) = [ (X8147X4 


- | (х5;-х2)4 * (YS-Y4) 


ғ 


These functions can be expressed as a column vector F(8) 
with transpose, 


Е' (0) = [ SR, 1+SR¡2+SR¡ 3+SRo] + -¿-+DSR]]+DSR, ),DSR, 31 +++] 


For m survey points F(0) has dimensions 3m + 3(m-1) = 3(2m-1). 


If F(0) could be measured without errors then some minimum 


BON: 





O u _ 


number of measurements would be sufficient to determine 9 
exactly. In this case the slant range differences are linear 
combinations of the slant ranges and the minimum data case, as 
shown by Vanderkulk (1961), requires 18 slant range measurements: 

SR; 4 (9), DE MD о 9 1.2.3 
or equivalently, 3 slant range measurements and 15 slant range 
differences 

SR, 5 (8). jee 23 

Deus ООо з= 1,2,3. 


If only slant range differences can be measured then 24 measure- 


ments, 


Il 


DSR. . (O9), i 
ANS 


ео, 29107253 
are sufficient to determine 0 exactly provided one has reasonable 
initial estimates of the parameters. This condition is a result 
of the linear approximation used to solve the non-linear equations. 

The observations Y differ from F(0) due to the presence 
of noise. Slant range observations from the transponder system 
can be expressed as 

Sij = 4. + SE 

where errors, Сад, аге assumed to ре independent identically 
distributed normal random variables with variance c Obser- 
vations of slant range differences from the Doppler system can 
be expressed as 


D$; - DSR; (9) + д: 


ШБ 





| 
| 


| 


where the errors, Ó is are assumed to be independent identi- 
cally distributed normal random variables with variance Ср. 
The errors for the complete system can be expressed as 


9 = [у - ғ(9)|. 


ШЕ Сһе еххогв Є апа 6 are independent of each other, than e 


10 


is characterized by the moment matrix 
Г Z 2 
0 “0427 
2 
"m Er 


I2 »3(т=1) 


5) Gauss-Newton Method of Non-linear Least Squares Estimation 

We have seen that the maximum likelihood estimator under 
conditions of additive, uncorrelated, Gaussian noise with zero 
mean is equivalent to the weighted least squares estimator and 
that this estimator is consistent and efficient. The weighted 
Ш ЕСЕ squares estimate is the 9 which minimizes the residual 
sum of squares 


en] Y [x-E0]. (7) 


— 


55(6) - e'V le - | 
Imn general the estimate Ə is the solution of the p normal 
equations 
Q SS (8) 
11-100 ep 
д Өз 


Eng that the solution is interior to the parameter space. 


When these equations are linear they can be solved directly. 
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When the normal equations are Ее as they are in this case 
an iterative technique must be used. A number of techniques 
exist; for example, quasi-Newton, Pen ace gradient, etc. A 
brief review of numerical techniques for fitting non-linear 
models is given by Chambers (1973). 

The method used to find the estimate Ө is based on a modi- 
fied Gauss-Newton procedure. The Gauss-Newton procedure consists 
ENunearly approximating the function F(8) and applying standard 
linear least squares techniques. 

One linear approximation of F(8) about an initial estimate 
8o is a simple Taylor series 


^ 


F(8) = Е(90) + x Ae 








dF 
where X = ES ^ is an n x p matrix of derivatives with 
= 
| ЖР» 
elements Xij = om and Ae = NE ES MSC O вата Vector 
A E 
un 90 


@eerength p. Choosing the new variable Z - Y - F (85) we have 


— 


the linear model 


-X A0 e. 


IN 


The linear approximation to the residual sum of squares 15 
^ 
ss(a) = evite = (2x Ae J "(zx Ael- ss(8y). 


The Gauss-Newton step Me is obtained Erzelimearsregress- 


ion techniques on this linear residual sum of' squares 
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(Draper and Smith, 1966); 


И 


nS) - 


Z 


As = (СМ ху” 
A new estimate 5, = бо +49 can be formed and the process 


repeated, until a given convergence criteria is reached. 


Mod tications of the Basic Gauss-Newton Procedure 


Modifications of the basic Gauss-Newton procedure can be 
made to insure that the step A98, will reduce the residual sum 
SE squares and avoid singularity problems in the х'м-іх 
matrix. A complete discussion of the modifications used, 
based on stepwise regression techniques, 1S contained in 
Jennrich and Sampson (1967). A summary of that discussion is 
presented here. 

One important modification of the basic Gauss-Newton 
 С(пге is the use of partial steps, NAS, in place of the 
full step A9 when the full step results in an increase in the 
residual sum of squares 55 (6). This frequently occurs when 
the linear approximations upon which the method is based fail. 
The proportion ^| to be used is obtained by trying the arbitrary 
Sequence of values 1,1/2, 1/4, 1/8,... until the residual sum 
Of squares is reduced. It can be shown that a sufficiently 
ШЕ! Әгер іп (һе direction of AS will always reduce the 


A 
residual sum of squares unless 9 is already a minimizing value. 
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Another modification of the basic procedure is dictated 
Meche possibility that the X'V^ix matrix is, or is nearly, 
singular. Step-wise regression techniques are used where the 
parameter selected at a given step is the one which makes the 
greatest reduction in the residual sum of squares. If there 
is a singularity problem then only a subset of the parameters 
is used. This also provides a convenient way to handle bound- 
ary restrictions on the parameter values. A parameter is 
modified only if the new parameter value is within or on the 


parameter boundary. 


wnt erative Technique 


In summary the procedure used consists of the following 
steps: 
the value of the function F(8) and its derivatives X 


are calculated using an initial estimate 85. 


1 


2) The X'V “X matrix is formed and the residual sum of 


squares SS (85) is computed. 


1 


3) The X'V X matrix is inverted in a step-wise manner 


1... -1. 1 


and the step size 46 = (X'V X) X'V ^Z is formed. 


4) A new estimate 91 = 85 + AO is formed and the new 
^ ^ 
residual sum of squares SS (81) computed. тө тте larger 
^^ р , ^ Z ^ 
Ш ЕП 55(90) then the step size is halved until SS (904) = 55 (87) 
ШЕ until Х = (1/2)+0. 


— 
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5) This completes an iteration. The process is repeated 
until the solution converges. -Convergence is determined by 


the relative change in the residual sum of squares, i.e. the 


solution has converged when the relative change is less than a 


specified criterion, C, 
A A 
S8(85,1) - SS(8,) 
^ 
S8 (9541) 

6) In practice we require that 4 successive iterations 
result in C'« C in order to accept the hypothesis of conver- 
gence. In some cases the algorithm may fail to converge to a 
minimizing solution (e.g. when it converges to a local, rather 
than a global, minimum). A plot of the parameter estimates 


and examination of the residuals will usually provide some 


Endscation of this. 


8) Error Estimates 


We have seen that maximum likelihood estimates are 
asymptotically normally distributed with known variances and 
that these variances are the diagonal.elements of the variance- 
covariance matrix I For a process perturbed by uncorrelated 
Gaussian noise with zero mean this matrix is given by the 


expression 


ЕО 





hec Bey tx) ee ш сыс Теа се correlation matrix and its 


inverse will be denoted by A. 


The standard error in estimating the ER parameter, oo 
1 


is the square root of the ¡th diagonal element of the variance- 
Seovariance matrix E 


Cg. = Oa; | (8) 
1 


th 


where a, ШЕ сһе i diagonal element of the A matrix 





n of, 2 = 
Wy. 


= à ^ | 
E Y 18 | 


| The terms on the right hand side of equation (8) are indepen- 
dent. The standard error, ©, depends upon the accuracy with 
nich the slant ranges and slant range ERE ees can be 
ШЕ... An estimate of the standard error can be obtained 
| from the residual sum of squares 

1 


Ga = ss (8) 
n-p 





where 55 (9) merderined rn equation (7). 

The other term, Va, 1S sometimes called the error 
| magnification factor and is a function of the survey geometry. 
Magnitude of the error magnification can be obtained by 


| evaluating equation (9). 
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The approximate errors in estimating the beacon coordin- 
ates 
A A А 
Jo. 8004 1 501002, 410 
can be reduced in two ways: 1) improvements in the measurement 
accuracy and 2) a suitable selection of the survey coordinates. 
In the next section, we will consider in detail how this may 


be accomplished. 
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III MINIMIZATION OF ERRORS 
1) Measurement errors ulse system 

Slant ranges and slant range differences are measured by 
two different systems for which the magnitude and sources of 
error will in general be different. A discussion of these 
errors requires a complete description of the system operations. 

The pulse system is designed to measure round trip travel 
time of an acoustic pulse between a ship-mounted or ship- 
suspended transducer and near-bottom transponders. When an 
acoustic pulse is transmitted from the transducer, digital 
counters in a shipboard receiver begin measuring elapsed time. 
The transponders receive and recognize the transmitted pulse 
and generate a "reply" pulse at а specific frequency. When 
a return pulse is detected by the shipboard receiver at a 
given frequency the corresponding counter is stopped and the 
elapsed time is displayed and transferred to a digital com- 
puter for processing. 

There are fixed time delays associated with this process; 
e.g. signal recognition time, Е turn around time 
Anar signal processing time. Accurate calibration of the system 


can reduce the uncertainties in each of these delay times to 


Mmeapproximately +.5 milliseconds (msec). For most applications 


paese timing errors can collectively be considered to be normally 
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distributed with zero mean and standard error of 1 msec. Two 
of these timing errors, Signal recognition time at the trans- 
ponder and at the ship are functions of the signal-to-noise 
ratio and the pulse bandwidth. For example, the errors in 
recognition time increase to approximately .8 to l msec for a 
signal-to-noise ratio at the pulse receiver of 27 dB. This 
signal-to-noise is typical of a range of 10 km in sea state 3 
for the system presently in use (source level: 189 dB, re: 

Up Pa pulse length: 10 msec). For ranges greater than about 

11 km (for the present system configuration) signal recognition 
is erratic and the errors can be assumed infinite. For fixed 
bandwidth, timing errors can be reduced by increasing signal 
power or they may be acounted for by some appropriate weighting 
of the data. 

Another source of error in the measurement of travel times 
is ship motion. The error in round trip travel time due to 
horizontal ship motion is given approximately by the expression 

d Mu MONA 


where = measured round trip travel time 


@ = angle from the vertical between the transponder and 
ship transducer 
Vx = horizontal component of ship speed away from the 


transponder 


Ac 





€ = average sound speed 
For ship speeds of 1 knot and measured travel times of 6 to 
12 seconds (typical in deep ocean applications) the round trip 
travel times will be in error by approximately 2 to 4 msec. 
This error can be reduced by taking measurements while stopped 
ог by reprocessing corrected data once’ the approximate ship- 
positions and speeds are known. Other ship motion such as 
heave, pitch and roll cannot easily be accounted for but when 
a large number of measurements are used their effect can be 
considered as an uncertainty in transducer depth equivalent to 
a timing error of approximately l msec. 

A summary of these errors is present in Table 1. To achieve 
these errors in measured travel time of a few milliseconds 
it is essential to use an accurately calibrated system and 
make the measurements while stopped. These errors can be 
assumed to have zero mean for a large number of measurements. 

To obtain geometric slant range from corrected travel 
time the velocity of sound must be known. Although the sound 
velocity has spatial and temporal EUER, horizontal and 
temporal variations can be assumed negligible. Sound velocity 
variation over depth, however, is commonly on the order of 
50 meters/sec for 5 km deep water. The structure of this 


variation, the sound velocity profile, can be determined from 


ааб 





Table 1. Sources of errors in measurement of Travel Time 


Source of Error 


¡Recognition Time 


Miscellaneous Timing Errors 
Horizontal Motion 2 KASE) 


Other Motion 
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historical data or direct measurements. Historical data can 
be expected to be in error by 2.5 m/sec or more, whereas 
present measurement techniques have an accuracy of about 
125 m/sec. For a precision navigation system the use of direct 
measurements is strongly preferred. When the sound velocity 
profile is known the travel time of an acoustic Signal between 
two points can be determined and an effective sound speed, Vor 
сап be calculated. 
Ve = SR/TT 

When the travel time is measured using a pulse navigation 
system the effective sound speed is not known since it depends 
upon the positions of the transponders and ship transducer. 
The slant range, however, can be determined approximately by 
using a calculated range R, equal to the corrected one-way 
travel time TTo multiplied by the arithmetic mean sound 


velocity c: 


The depth of the transducer ZS is usually known and errors of 
less than 50 meters in the estimate of transponder depth, 74, 


have a negligible effect on C. The arithmetic mean sound 
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Velocity can therefore be considered a on Constant. 

The ratio of slant range,. SR, to calculated range R = ет 
сап be expressed as a function K(Ve): 

K(Ve) = SR/R = Ve/c 
For a given sound velocity profile and depths ZS and 2; the 
effective sound speed, Ve, is a monotonically increasing 
function of the geometric angle between the transponder and 
the survey point measured from the vertical (Vass, 1966) (see 
Fig. 3). Since the calculated range, R, is also an increasing 
function of the geometric angle, K(Ve) can be expressed as a 
function of R and approximated by a second order polynomial 
K(Ve) = К(В) 2 Ө 2 

The coefficients, а), аге found by calculating the range В 
and the ratio SR/R for various angles and applying a second 
order curve fitting subroutine. Travel times are computed 
using standard ray acoustic techniques ar. T958) TA 
description of program ARCOR which performs these computations 
is contained in Appendix A. The error in using the relation- 
ship SR = R is usually less than "EC Бас can Бе ас 
great as 3 or 4 meters for angles greater than 80°. When the 
relationship SR = В (ар + а, Е + aR) is used the errors in 
Slant range estimates are less than 1 meters. When the sound 


velocity profile is accurately known, errors caused by acoustic 
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refraction are negligible compared to timing and transducer 


motion errors. 


2) Measurement Errors (Doppler system 
The Doppler system in its present configuration (Fig. 4) 
consists of three bottom-moored CW beacons each separated by 
50 Hz, a phase detector for each beacon and a real time com- 
MET processor. Each beacon contains a crystal oscillator 
with a frequency stability of 1073 Hz and has a nominal source 
level of 166 dB (re +1 Pascal @ 1 m). The beacon signal after 
being processed by its phase detector is recorded on analog mag- 
netic tape and is fed to a digital computer for real time pro- 
cessing. The phase detector output consists of two base-band 
outputs one proportional to the sine of the input phase angle, 
the other proportional to the cosine of the phase angle. Phase 
Variations due to a finite receiver velocity v(t), can be written 
EENS (t) = 2 Tr£ tv (t)/c. If the input signal is represented by 
s(t) = A(t) sup sli: = g(t) | | 
the phase detector outputs are given by 
51 (Е) = Att )pocos в (Е) | 
and s, (t) Ен И (Е) 5 
The instantaneous phase of the beacon signal can be determined 


From the ratio 


2,0 = tan Pl(t) - 


To- 
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practice the quadrant into which the phase angle falls is 
determined simply by examining the signs of 57 (6) апа s, (t). 
Phase differences between phase detector samples are computed 
and accumulated in quanta of А /4 Тос Ос ЕН”! The 
accumulated phase difference after T Seconds is of the form 
N | 

= ұм UN - gt. 1)| , N= Tf, 
where Е is the sampling frequency of the phase detector out- 
MS LO” Ez) (Porter et al., 1973). 

The distance, DSR, travelled toward or away from a beacon 
in time T is the accumulated phase Pr multiplied by the wave- 
length Aof the signal: 

DSR = 9, A 
The error in measurement of accumulated phase depends upon the 
number of samples N = Tf, and the Signal detection probabili- 
ties of the receiver. Since the samples can be considered to 
be independent measurements the variance of the accumulated 
phase is the sum of the variance of the individual samples 
2 N 2 2 
د‎ o mm 
n= n 


2 
and Cg = МР) 


ee 


th 


where Pe is the probability of error in the n sample. The 


probability of error is a function of the signal-to-noise | 
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ratio at the phase detector (see Appendix B). Table 2 shows 
some representative values of Бо їп а .1 sec interval for 
various signal-to-noise ratios encountered in practice. The 


associated errors Og and Og (T = 1 hr) are also given in 
AL 


n 
both number of quadrants and meters. A quarter wavelength of 
3 cm is assumed. Nominal horizontal range is calculated for 

a source level of 166 dB (re: l Pa) and a noise level of 

44 dB/Hz (re: 1 нра) assuming transmission losses are due only 
to spherical spreading (20 log r) and attenuation (1 dB/km). 
Sample calculations are given in Appendix C. 

Random phase fluctuations can be caused by multipath 
interference and forward scattering. The root mean square 
phase fluctuations from these sources has been estimated to 
EE out 13 quadrants (Porter et al., 1973). When compared 
with the errors due to ambient sea noise shown in Table 2 
they have little effect (about 5-159; increase іп Cg, ) and 
can be considered negligible. 

Slant range difference is ораса from accumulated phase 
change by multiplication by a scaling factor: the wavelength 
| Oof the signal. Errors in wavelength AX = c/f, are caused 
primarily by errors in sound e. instability ofthe 


beacon frequency is negligible). When the arithmetic mean 


Lon 








TABLE 2. Doppler Errors vs Signal-to-Noise Ratio 


10 log S/N | P | Og Ig (1 hr) 
Eu к ME | (meters)| 5! 
33 .023 | 217 41 | > | о 
| 253 | 48 а.а | 3.7 

мам ы ums 

| | | .293 Sie? 5% 
341 | 65 | 19 | 7.5. 
393 ЖИН 32252 
456 | 87 2.6 11. 
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sound velocity of an accurately measured sound velocity pro- 
ЕЕЕ 15 used, the standard error in wavelength will usually 
sS than 1 meter. This can be reduced further by using 
the average sound, 
5 = length of travel path/travel time, 

once the approximate survey positions are known. 

In summary, when operating above the pulse system thresh- 
E Slant ranges can be measured by the pulse system with an 
accuracy of 3 to 5 meters. Under the same conditions slant 
range differences can be measured with an accuracy of 1.2 to 


2.5 meters for measuring intervals of 1 hour. 


3) Selection of Survey Coordinates 


The error in estimating beacon coordinates also depends 

on the error magnification terms of the beacon coordinates, 
m= AAR Sal К 

where the a, are the oe diagonal elements of the A = x vix 
matrix. These terms are dependent upon the number and relative 
locations of the survey points. In general the error magnifi- 
cations will vary inversely with the degrees of freedom which 
increase with the number of survey points above a certain 
minimum. For example the slant range difference equations 
mequire at least 9 survey points for a solution to exist, and 


the degrees of freedom for this set of equations is equal to 


E 





E umber of additional survey points. Comparison of survey 
techniques should be based on the same number of degrees of 
Breedom. For example a pulse survey of 7 points (1 degree of 
freedom) should be compared to a Doppler survey of 10 points, 
and a pulse survey of 10 points should be compared to a pulse- 
Doppler survey of 10 points. 

A measure of "good" survey locations is the trace of the 
beacon coordinate covariance sub-matrix, that is, the sum of 
the squares of the error magnification factors of the beacon 


Coordinates: 


The trace is a function of the beacon and survey location 
coordinates, i.e. the parameter vector 0, and will be denoted 
by TR(9). For a given number of survey locations there will be 
a set or sets of survey locations for which the trace is a mini- 
Nn. A typical deep ocean deployment of the pulse-Doppler 
navigation system will be discussed for illustration. 

The most commonly used 3-beacon acoustic navigation net 
is ideally an equilateral triangle. The baseline length, i.e. 
mee length of the sides of the triangle, are typically on the 
order of the beacon depths: about 5 km. For this configura- 


tron of the net the covariance matrix was evaluated for several 


sets of survey locations, each consisting of 10 survey points. 
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A few examples are given in Figures 5A 2 5E. The value 
ST the trace of the beacon covariance matrix for each arrange- 
Meme 15 listed in Table 3. The trace for the corresponding 
survey with only pulse data is given for comparison. A more 
rigorous technique for determining "optimum" survey locations 
was also used. Powell's algorithm to find a (local) minimum 
of a function of several variables was applied to a function 
РИ) TR(O) + C(O), where C(9) is an arbitrary function added 
to the trace to constrain the survey locations within the 
maximum range of the system (Powell, 1965). No significant 
improvement to the starting value of G(8) © 3.6 was realized 


although several survey geometries were tried. 


4) Computer Simulation 


A computer program called PDSRV was developed to solve 
the survey problem using the parameter estimation techniques 
discussed. It is a modification of a Biomedical Computer Pro- 
gram (BMD 07R) (see Dixon, 1973). The program requires pulse 
data from at least six survey points and/or Doppler data from 
Beast ten survey points. In addition to the data, the 
user must input the average sound velocity and the acoustic 
merraction coefficients found by program ARCOR. The beacon 
frequencies, the depth(s) of the transducer and initial 


estimates of the beacon positions must also be provided. 
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TABLE 3. Comparison of the Trace of the Beacon Covariance Matrix 


for Pulse/Doppler and Pulse Surveys 


Figure # Trace Trace 
Pulse & Doppler Pulse Only 


5A 3.6 2776 


5D lO ы; ЕЕ 22 


5E 1968 225 
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Maximum and minimum values for the beacon coordinates are 
constrained to provide more rapid convergence. The constraints 


are set as follows: 


x. = X. +2000 m 
J J 

Y. = Y. +2000 m 
J J 

и 1000 mM 
J J 


where (х, У... 4.) are the initial estimates of the ¿th beacon 
coordinates. Other parameter bounds may be specified if desired. 
Although it is unlikely, even with the most rudimentary primary 
navigation system, that initial estimate errors will exceed 

those specified above. 

The program consists of three main sections: (1) initiali- 
zation and data input; (2) minimization of weighted squared 
error; (3) determination of the'approximate error covariance 
matrix of the beacon positions. For convenience the program 
is divided into 7 subroutines as illustrated in Figure 6. 

A description of each subroutine is given in Table 4. 

Two additional subroutines were used to simulate survey pro- 
blems. Subroutine READ calculates data for a given set of survey 
locations and beacon coordinates. Subroutine NORMA generates 


Gaussian distributed numbers with a specified mean and standard 


Stations. These are used to simulate data with errors. 
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' Name of 


| Subroutine 


= te — 


MAIN 


DINPT 


MINIZ 


XPRMX 


TABLE 4. 
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Description of Program PDSRV 


Description 


-eae ww - ч „= 


Reads the following user inputs: 


l1. Transducer depths 

Ze Convergenee erıterıa 

3. Bounds for beacon parameters 

4. Beacon coordinate estimates 

5. Survey position estimates (optional) 
6. Bounds of Survey position (optional) 


Reads input data from specified device(s) 
Calculates slant ranges (SR) from travel 
times 

Calculates changes in slant range (DSR) 
from accumulated phase change 


Computes the mean squared error (EE) 
from the weighted sum of squared errors: 
EE = 1/DOF) ss(8) 
Checks that EE is less than previous 
value; if not the step size is halved 
Checks convergence criteria 
-1 

Computes the elements of the X V X 
matrix 
Calculates the weighted sum of square 
errors: 

n 


ss(8) = Т7 (y.-£;) wi 
: 1=1 


256 
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TABLE 4. Continued 
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Subroutine | Description 
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| 1. Calculates values of the slant range 

| and slant range differences (£i) using 

| latest estimates of beacon and ship 

| . positions 
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1.е. the derivatives of the slant ranges 
and slant range differences: 
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STEP | fashion 
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Computes and outputs the following 
1. Final. estimates 
бота тон мае for Beacon Positions 
3. Root mean square error (standard erroe D 
ЛАШЕ етише сиса Ето terms of the 


DUTPT 
| covariance matrix 
| 5. Residuals and standard errors of each 
| measurement (opbtional) 
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Typical output of a simulated survey problem is shown in 
Tables 5a through 5d. Table 5a lists the beacon and survey 
coordinates used by subroutine READ and the errors in data 
generated by subroutine NORMA. Table 5b lists the simulated 
data which is input to program PDSRV. The final estimates, 
correlation matrix, error magnification factors and standard 
errors of the beacon coordinates are also shown in Table 5p. 
Similar information for the survey locations is listed as 
illustrated in Table 5c. The residuals, 

SR(observed) — SR(estimated) and 
DSR(observed) — DSR(estimated) 
can be compared with the estimates of standard error of each 
measurement. This output, labeled SURVEY ERRORS, is optional 
and is shown in Table 5c. 

After each iteration of the minimization process the mean 
Squared error and latest estimates of the beacon positions 
can be output as shown in Table 5d. The number of the iteration 
and the number of times the step size was halved is also 
included. 

Survey problems with measurement errors of more than 750 
meters have been simulated. The minimization process converged 
in every case and the beacon soria te estimates were nearly 


Bormally distributed about the true values with the standard 
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deviation approximately the н as that predicted. 

Single simulations with measurement Econo 
meters were made to determine the E uve. of the estimates 
of standard error, с, and the approximate error covariance 
matrix. Table 6 lists some examples showing the simulated 
measurement error, F the estimated standard error, o , and the 
estimate of the TRACE of the beacon error covariance matrix. 

The stability of the trace indicates the linearity of the model 
in the region around 9. 

The stability of the trace has another important implica- 
tion: suppose some "optimum" set of survey points have been 
determined for a particular navigation net; the actual survey 
points used might differ from the optimum by hundreds of meters 
due to navigation errors while conducting the survey; the 
resulting errors in the estimates of beacon coordinates, how- 
ever, would not be very different. The primary geometric factor 


EN Nffects the errors of the estimates is the distribution 


of the survey points not their exact location. 
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B DISCUSSION 

Using the parameter estimation technique and the procedures 
described, the relative positions of acoustic navigation beacons 
can be determined with an accuracy of Łl meter. 

A non-linear model is developed for a combined pulse-Doppler 
navigation system whose measurements of travel time and accumula- 
ted phase are perturbed by additive zero-mean Gaussian noise. 

The maximum likelihood estimation for the beacon coordinate: 

is found. The numerical technique used to solve the non-linear 
estimation problem is a modified Gauss-Newton method. This method 
consists of approximating the non-linear model with a linear 
expansion and iteratively solving the linear model using standard 
least Squares estimation techniques. An ee error co- 
variance matrix is found from which errors in the beacon coordi- 
nate can be estimated. Computer simulation has shown this 
technique to be stable for measurement errors within (and beyond) 
the range of errors encountered in practice. 

Many sources of error are present in the survey problem. 

The preceding sections detail the several steps to be taken in 
order to achieve an accuracy of +1 meter in the beacon estimates. 
These can be generalized as follows: 

im oound velocity should be а a precision 


instrument, for example, a Conductivity-Temperature-Depth (CTD) 


Bee 





mbe (see Gregg and Cox, 1971). Errors in the measurement of 
sound velocity will usually bias the estimates but will not 
affect their variance. 

2. Corrections should be applied to the measurements to 
E cunt for acoustic refraction. When these corrections are 
applied, the errors due to acoustic refraction can be consi- 
dered negligible. 

3. Travel time measurements should be made while the 
Surveying platform is dead-in-the-water, or nearly so, to 
eliminate errors due to horizontal platform motion. This is 
a potential source of large errors in the pulse system, although 
it does not affect the accuracy of the Doppler system measurement. 

4. The time interval between selected Doppler survey points 
should be less than 1 hour since errors in the measurement of 
accumulated phase increase with time. Intervals of approximately 
1/2 hour each are recommended. 

5. Survey locations which will minimize the error magnifi- 
cation factors of the beacon coordinate estimates should be 
used. The choice of survey locations is critical to obtaining 
accurate estimates and can also affect the ability of the 
estimation technique to converge to a global rather than a 
local minimum. Several examples of -— survey point selections 


Given in Section III. 3 above. 
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6. Pulse and Doppler measurements should be made simul- 
taneously at no less than 10 E cy locations. The use of 
more than 10 survey locations will, in general, reduce the 
errors in the beacon estimates. For example, the use of 40 
Eves locations can reduce these errors by a factor of 1/2. 
In most cases the capacity of the processing computer will limit 
the maximum number of survey locations which can be used 


Practıcally. 
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APPENDIX A 


DESCRIPTION OF PROGRAM ARCOR 

Pulse navigation systems calculate the distance between 
two points from the measured time of travel of an acoustic 
pulse between those two points. This distance can be approxi- 
mated by multiplying the travel time by an average sound 
velocity. This is an approximation since the sound travels 
in a refracted path rather than a straight line. The purpose 
БЕ Program ARCOR is to improve the approximation by accounting 
for the refraction of sound in water. For a given sound 
velocity profile and the depths of a sound source and receiver, 
this program uses standard ray-tracing techniques (Officer, 
1958) to find the time of travel, TT, of an acoustic pulse 
between the two points is the travel time TT multiplied by the 


arithmetic mean sound velocity c: 


The actual distance, SR is given by the expression 


2 107 2 


SR = (X| + 22) 
where X is the horizontal distance traveled and Z is the verti- 
cal distance between source and receiver. The ratio between 
the actual distance and the approximate distance, K; = SRi/Ri, 
is calculated for several travel times. The set of approximate 


distances and their ratios (В.,К; ), i = 1,N, are fit to a 
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second degree polynomial using standard least Squares re- 


gression: 


SR 2 
М ак = ак 
R 
where the a, are the acoustic refraction coefficients. The 


Metal distance SR can therefore be calculated when the 
approximate distance R = cTT is known: 
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О 1 2 

The set of approximate distances and their associated 
ratios are generated by calculating the necessary parameters 
for various angles 86 at which sound leaves the source (see 


Pig. 7). For each angle 99 there is a constant given by 


Snell's law 


= p (constant) = 
Co с; 
where су = sound velocity at the source and с. = sound velo- 


КОО at interface i. This constant and the sound velocity 


gradients 


determine the travel time and the distance an acoustic signal 
must travel to reach the receiver depth, assuming constant 
gradient in each layer. The angle at which sound leaves a 


layer 1s 


~ 
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ll 
9. = sin (pc). 


Шие travel time for each layer is 


1 | tan 0, /2 
Е = — іп (----). 
а tan 0: 1/2 
The horizontal distance traveled is 
1 
х = —  (cos0. йы ecco OS 95S) 
+ ^^ pg = 1 


For each starting angle 9 Che totalltravel time and horizontal 
distance traveled are found by summing the above quantities 
for each layer. We then have 


Total horizontal distance: X = D 


Vertical distance: 4 — Source depth-receiver depth 
Total travel time: TT = sit 

AN 1720 
Distance between source and receiver: SR= (X +42 ) 
Approximate distance: = етт 
Patio: K = SR/R 


The program computes the first starting angle by defining an 
angle, & = 90? and setting the departure angles equal to 

90° - X. Each successive “ is reduced by a factor, for 
example ,95, and a new departure angle is computed. The values 
of K and R are fitted with a second degree polynomial to pro- 
duce the coefficients, а. . A typical output is shown оп раде 


ES. A listing of the program begins on page 74. The pro- 


gram was originated by W.M. Marquet of the Woods Hole Oceanographic 
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Institution and includes calculations of coefficients for the 
case of a receiver at the same depth of the source. This 
Particular computation is not used in the survey problem and 
is not discussed here. The original version of this program 
was called Program SETUP and written by Mary Hunt and Roger 


Goldsmith of the Woods Hole Oceanographic Institution. 
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ARTIHMETIC MEAN SOUND VELOCITY = 15127304 
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Output of Program ARCOR. X=Horizontal distance between source and 
Bocojvor. SR = Slant range; TT — Travel timo; “ТЕТЛА Е лое о 
departure from source; к= SR/R; R= СТТ = approximate slant range. 
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APPENDIX B 


КООШО LON PROBABILITIES OF WOODS HOLE OCEANOGRAPHIC INSTITUTION 
NOTE COUNTER 

The detection probabilities of the Woods Hole Oceanographic 
Ion cycle counter can be calculated by considering the 
problem of detecting a sine wave in the presence of Gaussian 
noise. This treatment is sufficient because the estimate of 
phase is not based on the optimum estimator 





SCE) 
where sı (t) and 52 (Е) are the quadrature components of the 
Mco the phase detector; the estimate of phase is based 
on the signs of s} (t) and $2 (Е). The ее. probabilities 
of the sine component will be developed for illustration. 
The corresponding probabilities of the cosine component can 
be developed in a similar manner with identical results. 
Consider a sine wave in the absence of noise 
A A O E 
The amplitude A(t) varies slowly compared with the phase (t) 
ШІ сап be considered a constant. The phase can be assumed 


ШҒоспіу distributed 


p(8) = 1/27 ,|/e| < т 


сағ 





The probability density of the signal p(s) is therefore given 
by the expression: 
D 1 


MEN са» 


The phase detector output in the presence of noise will be 


S 





< f 








denoted by z(t) = s(t) + n(t) where the noise is normally 
distributed: 
2 
1 и О 
n no. | 
p(n) = © Qa, 7 r.m.s. noise voltage. 


TO, 


The probability density p(z) may be found by the convolu- 


tion of p(s) and p(n) assuming statistical independence of the 


two processes. This leads to the following integral for p(z): 
2 
1 -n^/20 1 1 
e EA 
“ZT OM | ғ 


This integral has several representations (see Ol'shevskii, 
1967) but is approximately a constant near z = 0 for A/J,>>1. 
Figure 8a illustrates the shape of this probability density 


2 
as a function of the parameter q? = 377 
n 


The probability of error is a weighted sum of the proba- 
bilities of false alarm P(FA) and false dismissal P(FD). Те 
cycle counter decides that the signal is positive when the 


voltage exceeds some threshold Vin > O.. Similarly the decision 


o 





that the signal is negative occurs when the voltage is less 
I ome negative threshold -Vp<0. A “false alarm" occurs 
when the signal is below the threshold and the signal plus 
noise exceeds the thresholds: 


|= (=) | > V_ given ls (2) | « V 


T T3? 


The probability of this event, P(FA), is approximately the same 
as the probability that the noise will exceed the thresholds 
in the absence of signal. This is represented by the shaded 


area in Figure 8b and is given by the expression 


РЕ ею 
P(FA) = 2 p(n)dn = 2 —— exp (-x*/2) dx. 
V Vin 5/27 2 
T Kan 


A false dismissal occurs when the signal is above the 
thresholds but the signal plus noise is below the thresholds: 
|| < Vp given Iso) > Up ANER roba bility OL this event 
P(FD), is approximately the same as the probability that the 
signal plus noise is below the thresholds for any signal. This 
is represented by the shaded area in Figure 8c and given 
by the expression 


Vip 
PFD) е f PAZ: 
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[mene region of interest p(z) can be approximated by a constant 
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Figure 8. Probability Densities of Signal Plus Noise 
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The probability of false dismissal can then be found 


NN = , 
ШШ ЧАК ЕУ г 
It can be seen that an increase in the threshold voltage- 
m noise ratio VO. will increase the probability of false 
dismissal and decrease the probability of false alarm for a 
EE M "onal-to-noise ratio. Since the losses assigned to 
these errors are approximately equal, the optimum threshold 


setting for a given signal-to-noise ratio will satisfy the 


following equation: 


P (FA) 
DS 
P (FD) 
or equivalently 2 5 
y 7X 72 УИС 
f 1T) e “ах------- 2. (10) 
V Тї (А/С. 
Lyon п) 


The expression on the left hand side of equation (10) is a 
tabulated function of Va C and can be solved by trial-and- 
error for various values of the signal power-to-noise ratio 
SN = ۸2/2 C. ^ (see Table 7). 

mié probability of error in determining the sign of the 
signal is given by the following expression for the optimum 
threshold setting: 
P, - 1/2 [P(FD) + P(FA)] = P(FD). (195) 
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10 log S/N 


TABLE 7 
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1.40 2150 
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EnNCerror in the cycle count, i.e. the number of quadrants 
of phase change between two successive cycles can be caused by 
an error in either or both of the quadrature components. When 
an error (false alarm or false dismissal) occurs in only one 
Esmponent the cycle count will be in error by +1 quadrant. 
when an error occurs An both components the cycle count will be 
Error by +2 quadrants. The following probability of cycle 
Seume error can be calculated from the probability of error ВЕ 
previously determined by equation (11): 


P(e-0) = (1-P,)? 


P(e=+l) = Елек | 


2 


Р(е=+2) = а , 


The variance in the nth sample can be calculated as follows 


= 2P 


2 
а № о. 


The variance of error in the measurement of accumulated phase 


change Øm after T seconds is 
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E 

| 

Y 
N 

> 
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Fr 


ОА 
d, | 


n 
where f, We the sampling frequency. Fora sampling frequency 
of 10 Hz and an interval of 1 hour = 3600 seconds the standard 


error in the measurement of Dep 15 © 251,90 Iy 
T n 


Eo 





Table 7 lists the values of Р, Ty and the optimum 
threshold-to-noise ratios for various signal-to-noise ratios 
typically encountered in practice. The variance is calculated 
assuming that the bias is zero or negligible. This assumption 


is valid when the velocity relative to a beacon is small. 
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DEPEND C 


SAMPLE CALCULATION OF SIGNAL-TO-NOISE RATIO 

In general the signal-to-noise level will decrease as the 
range between a sound source (e.g. CW beacon) and a receiver 
(e.g. ship transducer) increases. A sample calculation is 
made to demonstrate the realtionship between signal-to-noise 
E Xt the source (S/N), the range between source and receiver 
R, and signal-to-noise ratio at the receiver, (S/N)... The 
signal level of the CW beacons in the present configuration of 
Ec se-Doppler system is 166 dB (re: lyPa). A typical 
ambient sea noise level in sea state 3 is 44 dB/Hz (re: 1 м Ра) 
or 54 dB/10 Hz (see Urick, 1967). Тһе signal-to-noise level in 
a 10 Hz band is therefore 112 dB. а. losses are 
assumed due to spherical spreading and attenuation. The loss 
due to spherical spreading is 20 log R where R is the range in 
ШЕЕ. The attenuation of a 13 kHz signal is approximately 
1 dB/km. The signal-to-noise level at the receiver (S/N), iS 
therefore: 10 log(S/N), = 112 dB - 20 log R - R x 1077. For 
a horizontal range of 10 km and water depth of 5 km the range 


between a bottom-moored source and a surface receiver is 


approximately 


воз 





В = (52 + 192 17? x 10? SI x 10? m 


The (S/N) in decibels is therefore 


10 log (S/N) + Е: О 1.2 = 19.8 dB. 
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